1. Exploring L2A Reflectance¶

Summary

In this notebook we will open a netCDF4 file from the Earth Surface Minteral Dust Source Investigation (EMIT), specifically the Level 2A (L2A) Reflectance product. We will inspect the structure and plot the spectra of individual pixels and spatial coverage of a single scene. After that we will take advantage of the holoviews streams to build an interactive plot.

Background

The EMIT instrument is an imaging spectrometer that measures light in visible and infrared wavelengths. These measurements display unique spectral signatures that correspond to the composition on the Earth's surface. The EMIT mission focuses specifically on mapping the composition of minerals to better understand the effects of mineral dust throughout the Earth system and human populations now and in the future. More details about EMIT and its associated products can be found in the README.md and on the EMIT website.

The L2A Reflectance Product contains estimated surface reflectance. Surface reflectance is the fraction of incoming solar radiation reflected Earth's surface. Different materials reflect different proportions of radiation based upon their chemical composition, meaning that this information can be used to determine the composition of a target. In this guide you will learn how to plot a layer from the L2A reflectance spatially and look at the spectral curve associated with individual pixels, which can be used to identify targets.

Requirements

  • Set up Python Environment - See setup_instructions.md in the /setup/ folder
  • Download the required EMIT data - See setup_instructions.md in the /setup/ folder

Learning Objectives

  • How to open an EMIT .nc file as an xarray.Dataset
  • Apply the Geometry Lookup Table (GLT) to orthocorrect the image.
  • How to plot the spectra of pixels
  • How to plot specific bands as images

Tutorial Outline

1.1 Setup
1.2 Opening The Data
1.3 Plotting Data - Non-Orthorectified
1.4 Orthorectification
1.5 Plotting Data - Orthorectified 1.6 Saving Orthorectified Data
1.7 Interactive Spatial and Spectral Plots


1.1 Setup¶

Import the required Python libraries.

In [1]:
from osgeo import gdal
import numpy as np
import math
import xarray as xr
import geoviews as gv
import holoviews as hv
import hvplot.xarray
import netCDF4 as nc

Download the L2A Reflectance EMIT scene located here to your /data/ folder, then define an object representing the file path, like below.

In [2]:
fp = '../data/EMIT_L2A_RFL_001_20220903T163129_2224611_012.nc'

1.2 Opening EMIT Data¶

EMIT Data is distributed in a hierarchical NetCDF4 (.nc) format. Inside the netCDF4 file there are 3 groups, the root group containing reflectance values accross the downtrack, crosstrack, and bands dimensions, the sensor_band_parameters group containing the wavelength of each band center, and the full-width half maximum (FWHM) or bandwidth at half of the maximum amplitude, and the location group containing latitude and longitude values of each pixel as well as a geometric lookup table (GLT). The GLT is an orthorectified image that provides relative downtrack and crosstrack reference locations from the raw scene to facilitate fast projection of the dataset.

To access the .nc file we will use the xarray library. xarray only support non-hierarchical (flat) datasets, meaning that when loading a NetCDF into an xarray.Dataset, by default only the root group is added, the others have to be manually added. Since xarray does not recognize these other groups, the keys cannot be listed using that library. If we need to check what keys are present we can use the netcdf4 library. In the case of EMIT data, the reflectance group is the root group, and contains reflectance values and some metadata. As mentioned above, the other groups are sensor_band_parameters and location. Assuming we don't know the groups present, lets first explore the heirarchical structure of the EMIT data using the netcdf4 library.

In [3]:
ds = nc.Dataset(fp)
ds
Out[3]:
<class 'netCDF4._netCDF4.Dataset'>
root group (NETCDF4 data model, file format HDF5):
    ncei_template_version: NCEI_NetCDF_Swath_Template_v2.0
    summary: The Earth Surface Mineral Dust Source Investigation (EMIT) is an Earth Ventures-Instrument (EVI-4) Mission that maps the surface mineralogy of arid dust source regions via imaging spectroscopy in the visible and short-wave infrared (VSWIR). Installed on the International Space Station (ISS), the EMIT instrument is a Dyson imaging spectrometer that uses contiguous spectroscopic measurements from 410 to 2450 nm to resolve absoprtion features of iron oxides, clays, sulfates, carbonates, and other dust-forming minerals. During its one-year mission, EMIT will observe the sunlit Earth's dust source regions that occur within +/-52° latitude and produce maps of the source regions that can be used to improve forecasts of the role of mineral dust in the radiative forcing (warming or cooling) of the atmosphere.\n\nThis file contains L2A estimated surface reflectances         and geolocation data. Reflectance estimates are created using an Optimal Estimation technique - see ATBD for          details. Reflectance values are reported as fractions (relative to 1). 
    keywords: Imaging Spectroscopy, minerals, EMIT, dust, radiative forcing
    Conventions: CF-1.63
    sensor: EMIT (Earth Surface Mineral Dust Source Investigation)
    instrument: EMIT
    platform: ISS
    processing_version: V1.0
    institution: NASA Jet Propulsion Laboratory/California Institute of Technology
    license: https://science.nasa.gov/earth-science/earth-science-data/data-information-policy/
    naming_authority: LPDAAC
    date_created: 2022-11-14T09:50:54Z
    keywords_vocabulary: NASA Global Change Master Directory (GCMD) Science Keywords
    stdname_vocabulary: NetCDF Climate and Forecast (CF) Metadata Convention
    creator_name: Jet Propulsion Laboratory/California Institute of Technology
    creator_email: sarah.r.lundeen@jpl.nasa.gov
    creator_url: https://earth.jpl.nasa.gov/emit/
    project: Earth Surface Mineral Dust Source Investigation
    project_url: https://emit.jpl.nasa.gov/
    publisher_name: USGS LPDAAC
    publisher_url: https://lpdaac.usgs.gov
    publisher_email: lpdaac@usgs.gov
    identifier_product_doi_authority: https://doi.org
    flight_line: emit20220903t163129_o24611_s000
    time_coverage_start: 2022-09-03T16:31:29+0000
    time_coverage_end: 2022-09-03T16:31:41+0000
    software_build_version: 010603
    product_version: 01
    history: PGE Input files: radiance_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_rdn_b0106_v01.img, pixel_locations_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_loc_b0106_v01.img, observation_parameters_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_obs_b0106_v01.img, surface_model_config=/beegfs/store/emit/ops/repos/emit-sds-l2a/surface/surface_20221020.json
    easternmost_longitude: -62.5120945327963
    northernmost_latitude: -39.3067591475017
    westernmost_longitude: -61.236221412633064
    southernmost_latitude: -40.39610428069674
    spatialResolution: 0.000542232520256367
    spatial_ref: GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
    geotransform: [-6.25120945e+01  5.42232520e-04 -0.00000000e+00 -3.93067591e+01
 -0.00000000e+00 -5.42232520e-04]
    day_night_flag: Day
    title: EMIT L2A Surface Reflectance 60 m V001
    dimensions(sizes): downtrack(1280), crosstrack(1242), bands(285), ortho_y(2009), ortho_x(2353)
    variables(dimensions): float32 reflectance(downtrack, crosstrack, bands)
    groups: sensor_band_parameters, location

When inspecting the ds object, we can see the metadata, dimensions, variables, and groups. To specifically view the groups, append .groups.keys() to the object.

In [4]:
ds.groups.keys()
Out[4]:
dict_keys(['sensor_band_parameters', 'location'])

Now read the root reflectance group as an xarray.Dataset and preview it.

In [5]:
refl = xr.open_dataset(fp)
refl
Out[5]:
<xarray.Dataset>
Dimensions:      (downtrack: 1280, crosstrack: 1242, bands: 285)
Dimensions without coordinates: downtrack, crosstrack, bands
Data variables:
    reflectance  (downtrack, crosstrack, bands) float32 ...
Attributes: (12/38)
    ncei_template_version:             NCEI_NetCDF_Swath_Template_v2.0
    summary:                           The Earth Surface Mineral Dust Source ...
    keywords:                          Imaging Spectroscopy, minerals, EMIT, ...
    Conventions:                       CF-1.63
    sensor:                            EMIT (Earth Surface Mineral Dust Sourc...
    instrument:                        EMIT
    ...                                ...
    southernmost_latitude:             -40.39610428069674
    spatialResolution:                 0.000542232520256367
    spatial_ref:                       GEOGCS["WGS 84",DATUM["WGS_1984",SPHER...
    geotransform:                      [-6.25120945e+01  5.42232520e-04 -0.00...
    day_night_flag:                    Day
    title:                             EMIT L2A Surface Reflectance 60 m V001
xarray.Dataset
    • downtrack: 1280
    • crosstrack: 1242
    • bands: 285
      • reflectance
        (downtrack, crosstrack, bands)
        float32
        ...
        long_name :
        Surface Reflectance
        units :
        unitless
        [453081600 values with dtype=float32]
    • ncei_template_version :
      NCEI_NetCDF_Swath_Template_v2.0
      summary :
      The Earth Surface Mineral Dust Source Investigation (EMIT) is an Earth Ventures-Instrument (EVI-4) Mission that maps the surface mineralogy of arid dust source regions via imaging spectroscopy in the visible and short-wave infrared (VSWIR). Installed on the International Space Station (ISS), the EMIT instrument is a Dyson imaging spectrometer that uses contiguous spectroscopic measurements from 410 to 2450 nm to resolve absoprtion features of iron oxides, clays, sulfates, carbonates, and other dust-forming minerals. During its one-year mission, EMIT will observe the sunlit Earth's dust source regions that occur within +/-52° latitude and produce maps of the source regions that can be used to improve forecasts of the role of mineral dust in the radiative forcing (warming or cooling) of the atmosphere.\n\nThis file contains L2A estimated surface reflectances and geolocation data. Reflectance estimates are created using an Optimal Estimation technique - see ATBD for details. Reflectance values are reported as fractions (relative to 1).
      keywords :
      Imaging Spectroscopy, minerals, EMIT, dust, radiative forcing
      Conventions :
      CF-1.63
      sensor :
      EMIT (Earth Surface Mineral Dust Source Investigation)
      instrument :
      EMIT
      platform :
      ISS
      processing_version :
      V1.0
      institution :
      NASA Jet Propulsion Laboratory/California Institute of Technology
      license :
      https://science.nasa.gov/earth-science/earth-science-data/data-information-policy/
      naming_authority :
      LPDAAC
      date_created :
      2022-11-14T09:50:54Z
      keywords_vocabulary :
      NASA Global Change Master Directory (GCMD) Science Keywords
      stdname_vocabulary :
      NetCDF Climate and Forecast (CF) Metadata Convention
      creator_name :
      Jet Propulsion Laboratory/California Institute of Technology
      creator_email :
      sarah.r.lundeen@jpl.nasa.gov
      creator_url :
      https://earth.jpl.nasa.gov/emit/
      project :
      Earth Surface Mineral Dust Source Investigation
      project_url :
      https://emit.jpl.nasa.gov/
      publisher_name :
      USGS LPDAAC
      publisher_url :
      https://lpdaac.usgs.gov
      publisher_email :
      lpdaac@usgs.gov
      identifier_product_doi_authority :
      https://doi.org
      flight_line :
      emit20220903t163129_o24611_s000
      time_coverage_start :
      2022-09-03T16:31:29+0000
      time_coverage_end :
      2022-09-03T16:31:41+0000
      software_build_version :
      010603
      product_version :
      01
      history :
      PGE Input files: radiance_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_rdn_b0106_v01.img, pixel_locations_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_loc_b0106_v01.img, observation_parameters_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_obs_b0106_v01.img, surface_model_config=/beegfs/store/emit/ops/repos/emit-sds-l2a/surface/surface_20221020.json
      easternmost_longitude :
      -62.5120945327963
      northernmost_latitude :
      -39.3067591475017
      westernmost_longitude :
      -61.236221412633064
      southernmost_latitude :
      -40.39610428069674
      spatialResolution :
      0.000542232520256367
      spatial_ref :
      GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
      geotransform :
      [-6.25120945e+01 5.42232520e-04 -0.00000000e+00 -3.93067591e+01 -0.00000000e+00 -5.42232520e-04]
      day_night_flag :
      Day
      title :
      EMIT L2A Surface Reflectance 60 m V001

    We can see that the information read in only contains the root variable (reflectance) and attributes metadata, not those from the groups we previously listed. Using those group names we can read the other groups into their own xarray dataset object.

    Read in the sensor_band_parameters group as an xarray dataset and preview it.

    In [6]:
    wvl = xr.open_dataset(fp,group='sensor_band_parameters')
    wvl
    
    Out[6]:
    <xarray.Dataset>
    Dimensions:      (bands: 285)
    Dimensions without coordinates: bands
    Data variables:
        wavelengths  (bands) float32 381.0 388.4 395.8 ... 2.486e+03 2.493e+03
        fwhm         (bands) float32 8.415 8.415 8.415 8.415 ... 8.806 8.807 8.809
    xarray.Dataset
      • bands: 285
        • wavelengths
          (bands)
          float32
          ...
          long_name :
          Wavelength Centers
          units :
          nm
          array([ 381.00558,  388.4092 ,  395.81583, ..., 2478.153  , 2485.5386 ,
                 2492.9238 ], dtype=float32)
        • fwhm
          (bands)
          float32
          ...
          long_name :
          Full Width at Half Max
          units :
          nm
          array([8.415, 8.415, 8.415, ..., 8.806, 8.807, 8.809], dtype=float32)

      Now we can merge the two xarray datasets into a single dataset.

      In [7]:
      ds = xr.merge([refl,wvl])
      ds
      
      Out[7]:
      <xarray.Dataset>
      Dimensions:      (downtrack: 1280, crosstrack: 1242, bands: 285)
      Dimensions without coordinates: downtrack, crosstrack, bands
      Data variables:
          reflectance  (downtrack, crosstrack, bands) float32 ...
          wavelengths  (bands) float32 381.0 388.4 395.8 ... 2.486e+03 2.493e+03
          fwhm         (bands) float32 8.415 8.415 8.415 8.415 ... 8.806 8.807 8.809
      Attributes: (12/38)
          ncei_template_version:             NCEI_NetCDF_Swath_Template_v2.0
          summary:                           The Earth Surface Mineral Dust Source ...
          keywords:                          Imaging Spectroscopy, minerals, EMIT, ...
          Conventions:                       CF-1.63
          sensor:                            EMIT (Earth Surface Mineral Dust Sourc...
          instrument:                        EMIT
          ...                                ...
          southernmost_latitude:             -40.39610428069674
          spatialResolution:                 0.000542232520256367
          spatial_ref:                       GEOGCS["WGS 84",DATUM["WGS_1984",SPHER...
          geotransform:                      [-6.25120945e+01  5.42232520e-04 -0.00...
          day_night_flag:                    Day
          title:                             EMIT L2A Surface Reflectance 60 m V001
      xarray.Dataset
        • downtrack: 1280
        • crosstrack: 1242
        • bands: 285
          • reflectance
            (downtrack, crosstrack, bands)
            float32
            ...
            long_name :
            Surface Reflectance
            units :
            unitless
            [453081600 values with dtype=float32]
          • wavelengths
            (bands)
            float32
            381.0 388.4 ... 2.486e+03 2.493e+03
            long_name :
            Wavelength Centers
            units :
            nm
            array([ 381.00558,  388.4092 ,  395.81583, ..., 2478.153  , 2485.5386 ,
                   2492.9238 ], dtype=float32)
          • fwhm
            (bands)
            float32
            8.415 8.415 8.415 ... 8.807 8.809
            long_name :
            Full Width at Half Max
            units :
            nm
            array([8.415, 8.415, 8.415, ..., 8.806, 8.807, 8.809], dtype=float32)
        • ncei_template_version :
          NCEI_NetCDF_Swath_Template_v2.0
          summary :
          The Earth Surface Mineral Dust Source Investigation (EMIT) is an Earth Ventures-Instrument (EVI-4) Mission that maps the surface mineralogy of arid dust source regions via imaging spectroscopy in the visible and short-wave infrared (VSWIR). Installed on the International Space Station (ISS), the EMIT instrument is a Dyson imaging spectrometer that uses contiguous spectroscopic measurements from 410 to 2450 nm to resolve absoprtion features of iron oxides, clays, sulfates, carbonates, and other dust-forming minerals. During its one-year mission, EMIT will observe the sunlit Earth's dust source regions that occur within +/-52° latitude and produce maps of the source regions that can be used to improve forecasts of the role of mineral dust in the radiative forcing (warming or cooling) of the atmosphere.\n\nThis file contains L2A estimated surface reflectances and geolocation data. Reflectance estimates are created using an Optimal Estimation technique - see ATBD for details. Reflectance values are reported as fractions (relative to 1).
          keywords :
          Imaging Spectroscopy, minerals, EMIT, dust, radiative forcing
          Conventions :
          CF-1.63
          sensor :
          EMIT (Earth Surface Mineral Dust Source Investigation)
          instrument :
          EMIT
          platform :
          ISS
          processing_version :
          V1.0
          institution :
          NASA Jet Propulsion Laboratory/California Institute of Technology
          license :
          https://science.nasa.gov/earth-science/earth-science-data/data-information-policy/
          naming_authority :
          LPDAAC
          date_created :
          2022-11-14T09:50:54Z
          keywords_vocabulary :
          NASA Global Change Master Directory (GCMD) Science Keywords
          stdname_vocabulary :
          NetCDF Climate and Forecast (CF) Metadata Convention
          creator_name :
          Jet Propulsion Laboratory/California Institute of Technology
          creator_email :
          sarah.r.lundeen@jpl.nasa.gov
          creator_url :
          https://earth.jpl.nasa.gov/emit/
          project :
          Earth Surface Mineral Dust Source Investigation
          project_url :
          https://emit.jpl.nasa.gov/
          publisher_name :
          USGS LPDAAC
          publisher_url :
          https://lpdaac.usgs.gov
          publisher_email :
          lpdaac@usgs.gov
          identifier_product_doi_authority :
          https://doi.org
          flight_line :
          emit20220903t163129_o24611_s000
          time_coverage_start :
          2022-09-03T16:31:29+0000
          time_coverage_end :
          2022-09-03T16:31:41+0000
          software_build_version :
          010603
          product_version :
          01
          history :
          PGE Input files: radiance_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_rdn_b0106_v01.img, pixel_locations_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_loc_b0106_v01.img, observation_parameters_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_obs_b0106_v01.img, surface_model_config=/beegfs/store/emit/ops/repos/emit-sds-l2a/surface/surface_20221020.json
          easternmost_longitude :
          -62.5120945327963
          northernmost_latitude :
          -39.3067591475017
          westernmost_longitude :
          -61.236221412633064
          southernmost_latitude :
          -40.39610428069674
          spatialResolution :
          0.000542232520256367
          spatial_ref :
          GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
          geotransform :
          [-6.25120945e+01 5.42232520e-04 -0.00000000e+00 -3.93067591e+01 -0.00000000e+00 -5.42232520e-04]
          day_night_flag :
          Day
          title :
          EMIT L2A Surface Reflectance 60 m V001

        1.3 Plotting Data - Non-Orthorectified¶

        Pick a random downtrack and crosstrack location. Here we chose 660, 370 (downtrack,crosstrack). Next use the isel() function from xarray and the hvplot.line() functions to first select the spatial position and then plot a line showing the reflectance at that location.

        In [8]:
        ds['reflectance'].isel(downtrack=660,crosstrack=370).hvplot.line(y='reflectance',x='bands', color='black')
        
        Out[8]:

        We can see some flat regions in the spectral curve around bands 127 - 141 and 187 - 212. These are where water absoption features in these regions were removed. Typically this data is noisy due to the moisture present in the atmosphere; therefore, these spectral regions offer little information about targets and can be excluded from calculations.

        Although they have been reassigned a value of -0.01, we can mask them to improve visualization, by using the where() function to select regions of the dataset where the reflectance value is not equal to -0.01.

        In [9]:
        ds['reflectance'] = ds['reflectance'].where(ds['reflectance']!=-0.01)
        ds
        
        Out[9]:
        <xarray.Dataset>
        Dimensions:      (downtrack: 1280, crosstrack: 1242, bands: 285)
        Dimensions without coordinates: downtrack, crosstrack, bands
        Data variables:
            reflectance  (downtrack, crosstrack, bands) float32 0.02187 ... -0.004341
            wavelengths  (bands) float32 381.0 388.4 395.8 ... 2.486e+03 2.493e+03
            fwhm         (bands) float32 8.415 8.415 8.415 8.415 ... 8.806 8.807 8.809
        Attributes: (12/38)
            ncei_template_version:             NCEI_NetCDF_Swath_Template_v2.0
            summary:                           The Earth Surface Mineral Dust Source ...
            keywords:                          Imaging Spectroscopy, minerals, EMIT, ...
            Conventions:                       CF-1.63
            sensor:                            EMIT (Earth Surface Mineral Dust Sourc...
            instrument:                        EMIT
            ...                                ...
            southernmost_latitude:             -40.39610428069674
            spatialResolution:                 0.000542232520256367
            spatial_ref:                       GEOGCS["WGS 84",DATUM["WGS_1984",SPHER...
            geotransform:                      [-6.25120945e+01  5.42232520e-04 -0.00...
            day_night_flag:                    Day
            title:                             EMIT L2A Surface Reflectance 60 m V001
        xarray.Dataset
          • downtrack: 1280
          • crosstrack: 1242
          • bands: 285
            • reflectance
              (downtrack, crosstrack, bands)
              float32
              0.02187 0.0235 ... -0.004341
              long_name :
              Surface Reflectance
              units :
              unitless
              array([[[ 2.18727179e-02,  2.34996825e-02,  2.24097315e-02, ...,
                        5.09486720e-02,  4.14011404e-02,  3.83834094e-02],
                      [ 3.33263353e-02,  3.24866250e-02,  3.09991669e-02, ...,
                        4.79816645e-02,  3.49218473e-02,  3.35670225e-02],
                      [ 2.82608774e-02,  2.95085441e-02,  3.07566449e-02, ...,
                        4.54613604e-02,  3.96514796e-02,  3.74362692e-02],
                      ...,
                      [ 1.96052473e-02,  2.25020908e-02,  2.39406414e-02, ...,
                        3.26444837e-03,  3.16202128e-03, -7.51453126e-03],
                      [ 2.00065393e-02,  2.06522215e-02,  2.49008406e-02, ...,
                        1.20231677e-02,  3.51698720e-03, -8.72141402e-03],
                      [ 2.71530263e-02,  3.27588655e-02,  3.27786207e-02, ...,
                        9.27224942e-03, -1.76297978e-03, -1.69551857e-02]],
              
                     [[ 2.42982004e-02,  2.82031763e-02,  2.50691939e-02, ...,
                        5.81472628e-02,  4.51918021e-02,  4.17384505e-02],
                      [ 3.65232043e-02,  3.49093042e-02,  3.46683189e-02, ...,
                        5.50228432e-02,  4.55729216e-02,  3.83390896e-02],
                      [ 3.30173150e-02,  3.38553526e-02,  3.21070403e-02, ...,
                        4.87715304e-02,  3.87773365e-02,  3.41869704e-02],
              ...
                        7.07050832e-03, -6.94668619e-04, -5.65717695e-03],
                      [ 3.96752320e-02,  3.97593267e-02,  3.80864665e-02, ...,
                        6.65917993e-03,  2.68160063e-03, -3.92679824e-03],
                      [ 4.98440377e-02,  4.66296859e-02,  4.63311598e-02, ...,
                        5.61159104e-03,  6.36478537e-04, -4.51020850e-03]],
              
                     [[ 3.73413414e-02,  3.84780802e-02,  3.48522291e-02, ...,
                       -2.28073983e-03,  9.28868030e-05, -8.90796597e-04],
                      [ 5.10101579e-02,  5.17280884e-02,  4.85933349e-02, ...,
                        2.90575484e-03,  1.59259536e-03, -2.72730878e-03],
                      [ 5.10404482e-02,  4.94877398e-02,  4.80912179e-02, ...,
                        1.53037999e-03,  3.42500047e-03, -4.19431413e-03],
                      ...,
                      [ 4.42151539e-02,  4.10754755e-02,  4.14982289e-02, ...,
                        1.40570458e-02,  5.50026353e-03, -6.35674223e-04],
                      [ 4.03680094e-02,  3.81549597e-02,  4.02985513e-02, ...,
                        7.66669540e-03,  1.24720263e-03, -6.69119321e-03],
                      [ 4.97166365e-02,  4.74475846e-02,  4.77834642e-02, ...,
                        1.22473929e-02,  2.73120869e-03, -4.34138905e-03]]],
                    dtype=float32)
            • wavelengths
              (bands)
              float32
              381.0 388.4 ... 2.486e+03 2.493e+03
              long_name :
              Wavelength Centers
              units :
              nm
              array([ 381.00558,  388.4092 ,  395.81583, ..., 2478.153  , 2485.5386 ,
                     2492.9238 ], dtype=float32)
            • fwhm
              (bands)
              float32
              8.415 8.415 8.415 ... 8.807 8.809
              long_name :
              Full Width at Half Max
              units :
              nm
              array([8.415, 8.415, 8.415, ..., 8.806, 8.807, 8.809], dtype=float32)
          • ncei_template_version :
            NCEI_NetCDF_Swath_Template_v2.0
            summary :
            The Earth Surface Mineral Dust Source Investigation (EMIT) is an Earth Ventures-Instrument (EVI-4) Mission that maps the surface mineralogy of arid dust source regions via imaging spectroscopy in the visible and short-wave infrared (VSWIR). Installed on the International Space Station (ISS), the EMIT instrument is a Dyson imaging spectrometer that uses contiguous spectroscopic measurements from 410 to 2450 nm to resolve absoprtion features of iron oxides, clays, sulfates, carbonates, and other dust-forming minerals. During its one-year mission, EMIT will observe the sunlit Earth's dust source regions that occur within +/-52° latitude and produce maps of the source regions that can be used to improve forecasts of the role of mineral dust in the radiative forcing (warming or cooling) of the atmosphere.\n\nThis file contains L2A estimated surface reflectances and geolocation data. Reflectance estimates are created using an Optimal Estimation technique - see ATBD for details. Reflectance values are reported as fractions (relative to 1).
            keywords :
            Imaging Spectroscopy, minerals, EMIT, dust, radiative forcing
            Conventions :
            CF-1.63
            sensor :
            EMIT (Earth Surface Mineral Dust Source Investigation)
            instrument :
            EMIT
            platform :
            ISS
            processing_version :
            V1.0
            institution :
            NASA Jet Propulsion Laboratory/California Institute of Technology
            license :
            https://science.nasa.gov/earth-science/earth-science-data/data-information-policy/
            naming_authority :
            LPDAAC
            date_created :
            2022-11-14T09:50:54Z
            keywords_vocabulary :
            NASA Global Change Master Directory (GCMD) Science Keywords
            stdname_vocabulary :
            NetCDF Climate and Forecast (CF) Metadata Convention
            creator_name :
            Jet Propulsion Laboratory/California Institute of Technology
            creator_email :
            sarah.r.lundeen@jpl.nasa.gov
            creator_url :
            https://earth.jpl.nasa.gov/emit/
            project :
            Earth Surface Mineral Dust Source Investigation
            project_url :
            https://emit.jpl.nasa.gov/
            publisher_name :
            USGS LPDAAC
            publisher_url :
            https://lpdaac.usgs.gov
            publisher_email :
            lpdaac@usgs.gov
            identifier_product_doi_authority :
            https://doi.org
            flight_line :
            emit20220903t163129_o24611_s000
            time_coverage_start :
            2022-09-03T16:31:29+0000
            time_coverage_end :
            2022-09-03T16:31:41+0000
            software_build_version :
            010603
            product_version :
            01
            history :
            PGE Input files: radiance_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_rdn_b0106_v01.img, pixel_locations_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_loc_b0106_v01.img, observation_parameters_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_obs_b0106_v01.img, surface_model_config=/beegfs/store/emit/ops/repos/emit-sds-l2a/surface/surface_20221020.json
            easternmost_longitude :
            -62.5120945327963
            northernmost_latitude :
            -39.3067591475017
            westernmost_longitude :
            -61.236221412633064
            southernmost_latitude :
            -40.39610428069674
            spatialResolution :
            0.000542232520256367
            spatial_ref :
            GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
            geotransform :
            [-6.25120945e+01 5.42232520e-04 -0.00000000e+00 -3.93067591e+01 -0.00000000e+00 -5.42232520e-04]
            day_night_flag :
            Day
            title :
            EMIT L2A Surface Reflectance 60 m V001

          Since these datasets are large, we can go ahead and delete objects we won't be using to conserve memory.

          In [10]:
          del refl
          del wvl
          

          Plot the filtered reflectance values using the same downtrack and crosstrack position as above.

          In [11]:
          ds['reflectance'].isel(downtrack=660,crosstrack=370).hvplot.line(y='reflectance',x='wavelengths', color='black')
          
          Out[11]:

          Without the noisy data we can better interpret the spectral curve and hvplot will do a better job automatically scaling our axes.

          We can also plot the data spatially. Find the band nearest the 850nm wavelength in the NIR, then plot the data spatially using the isel() function to select only that band and using hvplot.image() to view the reflectance at 850nm of each pixel across the acquired region.

          In [12]:
          b850 = np.nanargmin(abs(ds['wavelengths'].values-850)) # Find band nearest to value of 650 nm (red)
          ds.isel(bands=b850).hvplot.image(cmap='viridis', aspect = 'equal', rasterize=True) 
          
          Out[12]:

          As previously mentioned a Geometry Lookup Table (GLT) is included in the location group of the netCDF4 file. Applying the GLT will orthorectify the image and give us Latitude and Longitude positional information.


          1.4 Orthorectification¶

          Before orthocorrecting, examine the location group from the dataset by reading it into xarray.

          In [13]:
          loc = xr.open_dataset(fp,group='location')
          loc
          
          Out[13]:
          <xarray.Dataset>
          Dimensions:  (downtrack: 1280, crosstrack: 1242, ortho_y: 2009, ortho_x: 2353)
          Dimensions without coordinates: downtrack, crosstrack, ortho_y, ortho_x
          Data variables:
              lat      (downtrack, crosstrack) float64 ...
              lon      (downtrack, crosstrack) float64 ...
              elev     (downtrack, crosstrack) float64 ...
              glt_x    (ortho_y, ortho_x) float64 ...
              glt_y    (ortho_y, ortho_x) float64 ...
          xarray.Dataset
            • downtrack: 1280
            • crosstrack: 1242
            • ortho_y: 2009
            • ortho_x: 2353
              • lat
                (downtrack, crosstrack)
                float64
                ...
                long_name :
                Longitude (WGS-84)
                units :
                degrees east
                [1589760 values with dtype=float64]
              • lon
                (downtrack, crosstrack)
                float64
                ...
                long_name :
                Latitude (WGS-84)
                units :
                degrees north
                [1589760 values with dtype=float64]
              • elev
                (downtrack, crosstrack)
                float64
                ...
                long_name :
                Surface Elevation
                units :
                m
                [1589760 values with dtype=float64]
              • glt_x
                (ortho_y, ortho_x)
                float64
                ...
                long_name :
                GLT Sample Lookup
                units :
                pixel location
                [4727177 values with dtype=float64]
              • glt_y
                (ortho_y, ortho_x)
                float64
                ...
                long_name :
                GLT Line Lookup
                units :
                pixel location
                [4727177 values with dtype=float64]

            We can see that each downtrack and crosstrack position has a latitude, longitude, and elevation, and the ortho_x and ortho_y data make up glt_x and glt_y arrays with a different shape. These arrays contain crosstrack and downtrack index values to quickly reproject the data. We will use these indexes to build an array of 2009x2353x285 (lat,lon,bands), filling it with the data from the EMIT dataset using the included emit_tools Python module. emit_tools contains some helpful functions for working with EMIT data using xarray.

            Import the emit_tools module and call use the help function to see how it can be used.

            Note: This function currently works with L1B Radiance and L2A Reflectance Data.

            In [14]:
            import sys
            sys.path.append('../modules/')
            from emit_tools import emit_xarray
            help(emit_xarray)
            
            Help on function emit_xarray in module emit_tools:
            
            emit_xarray(filepath, ortho=True, qmask=None, unpacked_bmask=None, GLT_NODATA_VALUE=0, fill_value=-9999)
                This function utilizes other functions in this module to streamline opening an EMIT dataset as an xarray.Dataset.
                
                Parameters:
                filepath: a filepath to an EMIT netCDF file
                ortho: True or False, whether to orthorectify the dataset or leave in crosstrack/downtrack coordinates.
                qmask: a numpy array output from the quality_mask function used to mask pixels based on quality flags selected in that function. Any non-orthorectified array with the proper crosstrack and downtrack dimensions can also be used.
                unpacked_bmask: a numpy array from  the band_mask function that can be used to mask band-specific pixels that have been interpolated.
                GLT_NODATA_VALUE: no data value for the GLT tables, 0 by default
                fill_value: the fill value for EMIT datasets, -9999 by default
                                
                Returns:
                out_xr: an xarray.Dataset constructed based on the parameters provided.
            
            

            We can see that the emit_xarray function will automatically apply the GLT to orthorectify the data unless ortho = False. The function will also apply masks if desired during construction of the output xarray.Dataset.

            Use the emit_xarray function to read in and orthorectify the L2A reflectance data.

            For a detailed walkthrough of the orthorectification process using the GLT see section 2 of the How_to_Orthorectify.ipynb in the how-tos folder.

            In [15]:
            ds_geo = emit_xarray(fp, ortho=True)
            ds_geo
            
            Out[15]:
            <xarray.Dataset>
            Dimensions:      (latitude: 2009, longitude: 2353, bands: 285)
            Coordinates:
              * latitude     (latitude) float64 -39.31 -39.31 -39.31 ... -40.39 -40.4 -40.4
              * longitude    (longitude) float64 -62.51 -62.51 -62.51 ... -61.24 -61.24
                wavelengths  (bands) float32 381.0 388.4 395.8 ... 2.486e+03 2.493e+03
                fwhm         (bands) float32 8.415 8.415 8.415 8.415 ... 8.806 8.807 8.809
                spatial_ref  int64 0
            Dimensions without coordinates: bands
            Data variables:
                reflectance  (latitude, longitude, bands) float32 nan nan nan ... nan nan
            Attributes: (12/38)
                ncei_template_version:             NCEI_NetCDF_Swath_Template_v2.0
                summary:                           The Earth Surface Mineral Dust Source ...
                keywords:                          Imaging Spectroscopy, minerals, EMIT, ...
                Conventions:                       CF-1.63
                sensor:                            EMIT (Earth Surface Mineral Dust Sourc...
                instrument:                        EMIT
                ...                                ...
                southernmost_latitude:             -40.39610428069674
                spatialResolution:                 0.000542232520256367
                spatial_ref:                       GEOGCS["WGS 84",DATUM["WGS_1984",SPHER...
                geotransform:                      [-6.25120945e+01  5.42232520e-04 -0.00...
                day_night_flag:                    Day
                title:                             EMIT L2A Surface Reflectance 60 m V001
            xarray.Dataset
              • latitude: 2009
              • longitude: 2353
              • bands: 285
              • latitude
                (latitude)
                float64
                -39.31 -39.31 ... -40.4 -40.4
                long_name :
                Latitude (WGS-84)
                units :
                degrees north
                array([-39.306759, -39.307301, -39.307844, ..., -40.394478, -40.39502 ,
                       -40.395562])
              • longitude
                (longitude)
                float64
                -62.51 -62.51 ... -61.24 -61.24
                long_name :
                Longitude (WGS-84)
                units :
                degrees east
                array([-62.512095, -62.511552, -62.51101 , ..., -61.237848, -61.237306,
                       -61.236764])
              • wavelengths
                (bands)
                float32
                381.0 388.4 ... 2.486e+03 2.493e+03
                long_name :
                Wavelength Centers
                units :
                nm
                array([ 381.00558,  388.4092 ,  395.81583,  403.2254 ,  410.638  ,
                        418.0536 ,  425.47214,  432.8927 ,  440.31726,  447.7428 ,
                        455.17035,  462.59888,  470.0304 ,  477.46292,  484.89743,
                        492.33292,  499.77142,  507.2099 ,  514.6504 ,  522.0909 ,
                        529.5333 ,  536.9768 ,  544.42126,  551.8667 ,  559.3142 ,
                        566.7616 ,  574.20905,  581.6585 ,  589.108  ,  596.55835,
                        604.0098 ,  611.4622 ,  618.9146 ,  626.36804,  633.8215 ,
                        641.2759 ,  648.7303 ,  656.1857 ,  663.6411 ,  671.09753,
                        678.5539 ,  686.0103 ,  693.4677 ,  700.9251 ,  708.38354,
                        715.84094,  723.2993 ,  730.7587 ,  738.2171 ,  745.6765 ,
                        753.1359 ,  760.5963 ,  768.0557 ,  775.5161 ,  782.97754,
                        790.4379 ,  797.89935,  805.36176,  812.8232 ,  820.2846 ,
                        827.746  ,  835.2074 ,  842.66986,  850.1313 ,  857.5937 ,
                        865.0551 ,  872.5176 ,  879.98004,  887.44147,  894.90393,
                        902.3664 ,  909.82886,  917.2913 ,  924.7538 ,  932.21625,
                        939.6788 ,  947.14026,  954.6027 ,  962.0643 ,  969.5268 ,
                        976.9883 ,  984.4498 ,  991.9114 ,  999.37286, 1006.8344 ,
                       1014.295  , 1021.7566 , 1029.2172 , 1036.6777 , 1044.1383 ,
                       1051.5989 , 1059.0596 , 1066.5201 , 1073.9797 , 1081.4404 ,
                       1088.9    , 1096.3597 , 1103.8184 , 1111.2781 , 1118.7368 ,
                ...
                       1796.4385 , 1803.8701 , 1811.3008 , 1818.7314 , 1826.1611 ,
                       1833.591  , 1841.0206 , 1848.4495 , 1855.8773 , 1863.3052 ,
                       1870.733  , 1878.16   , 1885.5869 , 1893.013  , 1900.439  ,
                       1907.864  , 1915.2892 , 1922.7133 , 1930.1375 , 1937.5607 ,
                       1944.9839 , 1952.4071 , 1959.8295 , 1967.2518 , 1974.6732 ,
                       1982.0946 , 1989.515  , 1996.9355 , 2004.355  , 2011.7745 ,
                       2019.1931 , 2026.6118 , 2034.0304 , 2041.4471 , 2048.865  ,
                       2056.2808 , 2063.6965 , 2071.1123 , 2078.5273 , 2085.9421 ,
                       2093.3562 , 2100.769  , 2108.1821 , 2115.5942 , 2123.0063 ,
                       2130.4175 , 2137.8289 , 2145.239  , 2152.6482 , 2160.0576 ,
                       2167.467  , 2174.8755 , 2182.283  , 2189.6904 , 2197.097  ,
                       2204.5034 , 2211.9092 , 2219.3147 , 2226.7195 , 2234.1233 ,
                       2241.5269 , 2248.9297 , 2256.3328 , 2263.7346 , 2271.1365 ,
                       2278.5376 , 2285.9387 , 2293.3386 , 2300.7378 , 2308.136  ,
                       2315.5342 , 2322.9326 , 2330.3298 , 2337.7263 , 2345.1216 ,
                       2352.517  , 2359.9126 , 2367.3071 , 2374.7007 , 2382.0935 ,
                       2389.486  , 2396.878  , 2404.2695 , 2411.6604 , 2419.0513 ,
                       2426.4402 , 2433.8303 , 2441.2183 , 2448.6064 , 2455.9944 ,
                       2463.3816 , 2470.7678 , 2478.153  , 2485.5386 , 2492.9238 ],
                      dtype=float32)
              • fwhm
                (bands)
                float32
                8.415 8.415 8.415 ... 8.807 8.809
                long_name :
                Full Width at Half Max
                units :
                nm
                array([8.415, 8.415, 8.415, 8.415, 8.417, 8.418, 8.419, 8.421, 8.422,
                       8.424, 8.425, 8.426, 8.428, 8.429, 8.431, 8.432, 8.433, 8.435,
                       8.436, 8.438, 8.439, 8.44 , 8.442, 8.443, 8.445, 8.446, 8.447,
                       8.449, 8.45 , 8.452, 8.453, 8.454, 8.456, 8.457, 8.459, 8.46 ,
                       8.461, 8.463, 8.464, 8.466, 8.467, 8.468, 8.47 , 8.471, 8.473,
                       8.474, 8.475, 8.477, 8.478, 8.48 , 8.481, 8.482, 8.484, 8.485,
                       8.487, 8.488, 8.489, 8.491, 8.492, 8.494, 8.495, 8.496, 8.498,
                       8.499, 8.501, 8.502, 8.503, 8.505, 8.506, 8.508, 8.509, 8.51 ,
                       8.512, 8.513, 8.515, 8.516, 8.517, 8.519, 8.52 , 8.522, 8.523,
                       8.524, 8.526, 8.527, 8.529, 8.53 , 8.531, 8.533, 8.534, 8.536,
                       8.537, 8.538, 8.54 , 8.541, 8.543, 8.544, 8.545, 8.547, 8.548,
                       8.55 , 8.551, 8.552, 8.554, 8.555, 8.557, 8.558, 8.559, 8.561,
                       8.562, 8.564, 8.565, 8.566, 8.568, 8.569, 8.571, 8.572, 8.573,
                       8.575, 8.576, 8.578, 8.579, 8.58 , 8.582, 8.583, 8.585, 8.586,
                       8.587, 8.589, 8.59 , 8.592, 8.593, 8.594, 8.596, 8.597, 8.599,
                       8.6  , 8.601, 8.603, 8.604, 8.606, 8.607, 8.608, 8.61 , 8.611,
                       8.613, 8.614, 8.615, 8.617, 8.618, 8.62 , 8.621, 8.622, 8.624,
                       8.625, 8.627, 8.628, 8.629, 8.631, 8.632, 8.634, 8.635, 8.636,
                       8.638, 8.639, 8.641, 8.642, 8.643, 8.645, 8.646, 8.648, 8.649,
                       8.65 , 8.652, 8.653, 8.655, 8.656, 8.657, 8.659, 8.66 , 8.662,
                       8.663, 8.664, 8.666, 8.667, 8.669, 8.67 , 8.671, 8.673, 8.674,
                       8.676, 8.677, 8.678, 8.68 , 8.681, 8.683, 8.684, 8.685, 8.687,
                       8.688, 8.69 , 8.691, 8.692, 8.694, 8.695, 8.697, 8.698, 8.699,
                       8.701, 8.702, 8.704, 8.705, 8.706, 8.708, 8.709, 8.711, 8.712,
                       8.714, 8.715, 8.716, 8.718, 8.719, 8.721, 8.722, 8.723, 8.725,
                       8.726, 8.727, 8.729, 8.73 , 8.732, 8.733, 8.734, 8.736, 8.737,
                       8.739, 8.74 , 8.741, 8.743, 8.744, 8.746, 8.747, 8.748, 8.75 ,
                       8.751, 8.753, 8.754, 8.755, 8.757, 8.758, 8.76 , 8.761, 8.763,
                       8.764, 8.765, 8.767, 8.768, 8.77 , 8.771, 8.772, 8.774, 8.775,
                       8.777, 8.778, 8.779, 8.781, 8.782, 8.783, 8.785, 8.786, 8.788,
                       8.789, 8.79 , 8.792, 8.793, 8.795, 8.796, 8.797, 8.799, 8.8  ,
                       8.802, 8.803, 8.804, 8.806, 8.807, 8.809], dtype=float32)
              • spatial_ref
                ()
                int64
                0
                crs_wkt :
                GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
                semi_major_axis :
                6378137.0
                semi_minor_axis :
                6356752.314245179
                inverse_flattening :
                298.257223563
                reference_ellipsoid_name :
                WGS 84
                longitude_of_prime_meridian :
                0.0
                prime_meridian_name :
                Greenwich
                geographic_crs_name :
                WGS 84
                grid_mapping_name :
                latitude_longitude
                spatial_ref :
                GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
                array(0)
              • reflectance
                (latitude, longitude, bands)
                float32
                nan nan nan nan ... nan nan nan nan
                long_name :
                Surface Reflectance
                units :
                unitless
                array([[[nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        ...,
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan]],
                
                       [[nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        ...,
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan]],
                
                       [[nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        ...,
                ...
                        ...,
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan]],
                
                       [[nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        ...,
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan]],
                
                       [[nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        ...,
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan],
                        [nan, nan, nan, ..., nan, nan, nan]]], dtype=float32)
            • ncei_template_version :
              NCEI_NetCDF_Swath_Template_v2.0
              summary :
              The Earth Surface Mineral Dust Source Investigation (EMIT) is an Earth Ventures-Instrument (EVI-4) Mission that maps the surface mineralogy of arid dust source regions via imaging spectroscopy in the visible and short-wave infrared (VSWIR). Installed on the International Space Station (ISS), the EMIT instrument is a Dyson imaging spectrometer that uses contiguous spectroscopic measurements from 410 to 2450 nm to resolve absoprtion features of iron oxides, clays, sulfates, carbonates, and other dust-forming minerals. During its one-year mission, EMIT will observe the sunlit Earth's dust source regions that occur within +/-52° latitude and produce maps of the source regions that can be used to improve forecasts of the role of mineral dust in the radiative forcing (warming or cooling) of the atmosphere.\n\nThis file contains L2A estimated surface reflectances and geolocation data. Reflectance estimates are created using an Optimal Estimation technique - see ATBD for details. Reflectance values are reported as fractions (relative to 1).
              keywords :
              Imaging Spectroscopy, minerals, EMIT, dust, radiative forcing
              Conventions :
              CF-1.63
              sensor :
              EMIT (Earth Surface Mineral Dust Source Investigation)
              instrument :
              EMIT
              platform :
              ISS
              processing_version :
              V1.0
              institution :
              NASA Jet Propulsion Laboratory/California Institute of Technology
              license :
              https://science.nasa.gov/earth-science/earth-science-data/data-information-policy/
              naming_authority :
              LPDAAC
              date_created :
              2022-11-14T09:50:54Z
              keywords_vocabulary :
              NASA Global Change Master Directory (GCMD) Science Keywords
              stdname_vocabulary :
              NetCDF Climate and Forecast (CF) Metadata Convention
              creator_name :
              Jet Propulsion Laboratory/California Institute of Technology
              creator_email :
              sarah.r.lundeen@jpl.nasa.gov
              creator_url :
              https://earth.jpl.nasa.gov/emit/
              project :
              Earth Surface Mineral Dust Source Investigation
              project_url :
              https://emit.jpl.nasa.gov/
              publisher_name :
              USGS LPDAAC
              publisher_url :
              https://lpdaac.usgs.gov
              publisher_email :
              lpdaac@usgs.gov
              identifier_product_doi_authority :
              https://doi.org
              flight_line :
              emit20220903t163129_o24611_s000
              time_coverage_start :
              2022-09-03T16:31:29+0000
              time_coverage_end :
              2022-09-03T16:31:41+0000
              software_build_version :
              010603
              product_version :
              01
              history :
              PGE Input files: radiance_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_rdn_b0106_v01.img, pixel_locations_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_loc_b0106_v01.img, observation_parameters_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_obs_b0106_v01.img, surface_model_config=/beegfs/store/emit/ops/repos/emit-sds-l2a/surface/surface_20221020.json
              easternmost_longitude :
              -62.5120945327963
              northernmost_latitude :
              -39.3067591475017
              westernmost_longitude :
              -61.236221412633064
              southernmost_latitude :
              -40.39610428069674
              spatialResolution :
              0.000542232520256367
              spatial_ref :
              GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
              geotransform :
              [-6.25120945e+01 5.42232520e-04 -0.00000000e+00 -3.93067591e+01 -0.00000000e+00 -5.42232520e-04]
              day_night_flag :
              Day
              title :
              EMIT L2A Surface Reflectance 60 m V001

            1.5 Plotting Data - Orthorectified¶

            Now that the data has been orthorectified, plot the georeferenced dataset using the same single wavelength (850nm) as above. We can use the aspect = 'equal' option to preserve the square pixel dimensions. The rasterize = True will help save memory and reduces the size of this notebook. For higher quality outputs, this can be omitted.

            In [16]:
            ds_geo.isel(bands=b850).hvplot.image(cmap='viridis', frame_width=500, aspect = 'equal', rasterize=True)
            
            Out[16]:

            We an also plot the data against an imagery tile using the geo=True and tiles= parameters instead of . Any tile source available in geoviews should work here. This will change the axis names, but that can be fixed by adding them manually in the options, like below.

            In [17]:
            ds_geo.isel(bands=b850).hvplot.image(cmap='viridis', frame_width=500, geo=True, tiles='EsriImagery',rasterize=True).opts(
                xlabel=f'{ds_geo.longitude.long_name} ({ds_geo.longitude.units})', ylabel=f'{ds_geo.latitude.long_name} ({ds_geo.latitude.units})')
            
            Out[17]:

            We can see that the orthorectification step rotated the image and placed it on a Lat/Lon grid. Now that we have a better idea of what the target area looks like, we can also plot the spectra using the georeferenced data. First, filter out the water absorption bands using the good_refl() function we wrote in section 1.3.

            In [18]:
            ds_geo = ds_geo.where(ds_geo['reflectance']!=-0.01)
            

            Now, plot the spectra at the Lat/Lon coordinates provided below.

            In [19]:
            ds_geo.sel(longitude=-61.833,latitude=-39.710,method='nearest').hvplot.line(y='reflectance',x='wavelengths', color='black', frame_width=400)
            
            Out[19]:

            1.6 Writing an Orthorectified Output¶

            At this point, the ds_geo orthorectified EMIT data can also be written as a flattened netCDF4 output that can be read using the xarray.open_dataset function, if desired.

            In [ ]:
            ds_geo.to_netcdf('../data/geo_ds_out.nc')
            # Example for Opening 
            # ds = xr.open_dataset('../data/geo_ds_out.nc')
            

            1.7 Interactive Spatial and Spectral Plots¶

            Combining the Spatial and Spectral information into a single visualization can be a powerful tool for exploring and inspecting data quality. Using the Streams function of Holoviews we can link a spatial map to a plot of spectra.

            We could plot a single band image as we previously have, but using an RGB image may help infer what targets we're examining. Build an RGB image following the steps below.

            Select bands to represent red, green, and blue by finding the nearest to a chosen wavelength.

            In [20]:
            # Find Nearest Bands
            b650 = np.nanargmin(abs(ds_geo['wavelengths'].values-650)) # Find band nearest to value of 650 nm (red)
            b560 = np.nanargmin(abs(ds_geo['wavelengths'].values-560)) # Find band nearest to value of 560 nm (green)
            b470 = np.nanargmin(abs(ds_geo['wavelengths'].values-470)) # Find band nearest to value of 470 nm (blue)
            

            Next, write a function to build an array from the chosen bands and scale the values using a gamma correction. Without applying this scaling the majority of the image would be very dark, with the reflectance data being skewed by the few pixels with very high reflectance.

            Note: This has no impact on analysis or data, just visualizing the RGB map.

            In [21]:
            def gamma_adjust(ds,band):
                # Define Array
                array = ds['reflectance'].isel(bands=band).values
                # Rescale Values using gamma to adjust brightness
                gamma = math.log(0.2)/math.log(np.nanmean(array)) # Create exponent for gamma scaling - can be adjusted by changing 0.2 
                scaled = np.power(array,gamma).clip(0,1) # Apply scaling and clip to 0-1 range
                scaled = np.nan_to_num(scaled, nan = 1) #Assign NA's to 1 so they appear white in plots
                return scaled
            

            Now apply this function to each of the selected bands, stack them, build the arrays of coordinates (Lat, Lon, Bands) needed to create an xarray.Dataset, then build the dataset.

            In [22]:
            # Scale the Bands
            r = gamma_adjust(ds_geo,b650)
            g = gamma_adjust(ds_geo,b560)
            b = gamma_adjust(ds_geo,b470)
            # Stack Bands and make an index
            rgb = np.stack([r,g,b]) # Stack r,g,b arrays and assign NA's to 1 so they appear white in plots
            bds = np.array([0,1,2])
            # Pull lat and lon values from geocorrected arrays
            x = ds_geo['longitude'].values
            y = ds_geo['latitude'].values
            # Create new rgb xarray data array.
            data_vars = {'RGB':(['bands','latitude','longitude'], rgb)} 
            coords = {'bands':(['bands'],bds), 'latitude':(['latitude'],y), 'longitude':(['longitude'],x)}
            attrs = ds_geo.attrs
            ds_rgb = xr.Dataset(data_vars=data_vars, coords=coords, attrs=attrs)
            ds_rgb.coords['latitude'].attrs = ds_geo['longitude'].attrs
            ds_rgb.coords['longitude'].attrs = ds_geo['latitude'].attrs
            ds_rgb
            
            Out[22]:
            <xarray.Dataset>
            Dimensions:    (bands: 3, latitude: 2009, longitude: 2353)
            Coordinates:
              * bands      (bands) int64 0 1 2
              * latitude   (latitude) float64 -39.31 -39.31 -39.31 ... -40.39 -40.4 -40.4
              * longitude  (longitude) float64 -62.51 -62.51 -62.51 ... -61.24 -61.24 -61.24
            Data variables:
                RGB        (bands, latitude, longitude) float32 1.0 1.0 1.0 ... 1.0 1.0 1.0
            Attributes: (12/38)
                ncei_template_version:             NCEI_NetCDF_Swath_Template_v2.0
                summary:                           The Earth Surface Mineral Dust Source ...
                keywords:                          Imaging Spectroscopy, minerals, EMIT, ...
                Conventions:                       CF-1.63
                sensor:                            EMIT (Earth Surface Mineral Dust Sourc...
                instrument:                        EMIT
                ...                                ...
                southernmost_latitude:             -40.39610428069674
                spatialResolution:                 0.000542232520256367
                spatial_ref:                       GEOGCS["WGS 84",DATUM["WGS_1984",SPHER...
                geotransform:                      [-6.25120945e+01  5.42232520e-04 -0.00...
                day_night_flag:                    Day
                title:                             EMIT L2A Surface Reflectance 60 m V001
            xarray.Dataset
              • bands: 3
              • latitude: 2009
              • longitude: 2353
              • bands
                (bands)
                int64
                0 1 2
                array([0, 1, 2])
              • latitude
                (latitude)
                float64
                -39.31 -39.31 ... -40.4 -40.4
                long_name :
                Longitude (WGS-84)
                units :
                degrees east
                array([-39.306759, -39.307301, -39.307844, ..., -40.394478, -40.39502 ,
                       -40.395562])
              • longitude
                (longitude)
                float64
                -62.51 -62.51 ... -61.24 -61.24
                long_name :
                Latitude (WGS-84)
                units :
                degrees north
                array([-62.512095, -62.511552, -62.51101 , ..., -61.237848, -61.237306,
                       -61.236764])
              • RGB
                (bands, latitude, longitude)
                float32
                1.0 1.0 1.0 1.0 ... 1.0 1.0 1.0 1.0
                array([[[1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        ...,
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.]],
                
                       [[1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        ...,
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.]],
                
                       [[1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        ...,
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.],
                        [1., 1., 1., ..., 1., 1., 1.]]], dtype=float32)
            • ncei_template_version :
              NCEI_NetCDF_Swath_Template_v2.0
              summary :
              The Earth Surface Mineral Dust Source Investigation (EMIT) is an Earth Ventures-Instrument (EVI-4) Mission that maps the surface mineralogy of arid dust source regions via imaging spectroscopy in the visible and short-wave infrared (VSWIR). Installed on the International Space Station (ISS), the EMIT instrument is a Dyson imaging spectrometer that uses contiguous spectroscopic measurements from 410 to 2450 nm to resolve absoprtion features of iron oxides, clays, sulfates, carbonates, and other dust-forming minerals. During its one-year mission, EMIT will observe the sunlit Earth's dust source regions that occur within +/-52° latitude and produce maps of the source regions that can be used to improve forecasts of the role of mineral dust in the radiative forcing (warming or cooling) of the atmosphere.\n\nThis file contains L2A estimated surface reflectances and geolocation data. Reflectance estimates are created using an Optimal Estimation technique - see ATBD for details. Reflectance values are reported as fractions (relative to 1).
              keywords :
              Imaging Spectroscopy, minerals, EMIT, dust, radiative forcing
              Conventions :
              CF-1.63
              sensor :
              EMIT (Earth Surface Mineral Dust Source Investigation)
              instrument :
              EMIT
              platform :
              ISS
              processing_version :
              V1.0
              institution :
              NASA Jet Propulsion Laboratory/California Institute of Technology
              license :
              https://science.nasa.gov/earth-science/earth-science-data/data-information-policy/
              naming_authority :
              LPDAAC
              date_created :
              2022-11-14T09:50:54Z
              keywords_vocabulary :
              NASA Global Change Master Directory (GCMD) Science Keywords
              stdname_vocabulary :
              NetCDF Climate and Forecast (CF) Metadata Convention
              creator_name :
              Jet Propulsion Laboratory/California Institute of Technology
              creator_email :
              sarah.r.lundeen@jpl.nasa.gov
              creator_url :
              https://earth.jpl.nasa.gov/emit/
              project :
              Earth Surface Mineral Dust Source Investigation
              project_url :
              https://emit.jpl.nasa.gov/
              publisher_name :
              USGS LPDAAC
              publisher_url :
              https://lpdaac.usgs.gov
              publisher_email :
              lpdaac@usgs.gov
              identifier_product_doi_authority :
              https://doi.org
              flight_line :
              emit20220903t163129_o24611_s000
              time_coverage_start :
              2022-09-03T16:31:29+0000
              time_coverage_end :
              2022-09-03T16:31:41+0000
              software_build_version :
              010603
              product_version :
              01
              history :
              PGE Input files: radiance_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_rdn_b0106_v01.img, pixel_locations_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_loc_b0106_v01.img, observation_parameters_file=/beegfs/store/emit/ops/data/acquisitions/20220903/emit20220903t163129/l1b/emit20220903t163129_o24611_s000_l1b_obs_b0106_v01.img, surface_model_config=/beegfs/store/emit/ops/repos/emit-sds-l2a/surface/surface_20221020.json
              easternmost_longitude :
              -62.5120945327963
              northernmost_latitude :
              -39.3067591475017
              westernmost_longitude :
              -61.236221412633064
              southernmost_latitude :
              -40.39610428069674
              spatialResolution :
              0.000542232520256367
              spatial_ref :
              GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]
              geotransform :
              [-6.25120945e+01 5.42232520e-04 -0.00000000e+00 -3.93067591e+01 -0.00000000e+00 -5.42232520e-04]
              day_night_flag :
              Day
              title :
              EMIT L2A Surface Reflectance 60 m V001

            Lastly, use the RGB data array to build a map object with the hvplot.rgb() function from holoviews.

            In [23]:
            # Define RGB Image
            map = ds_rgb.hvplot.rgb(x='longitude', y='latitude', bands='bands', aspect = 'equal', frame_width=400)
            

            To visualize the spectral and spatial data side-by-side, we use the pointerXY and Tap features of the streams functionality from the holoviews library. First, define objects resulting from the stream of the pointer x and y position on a spatial plot, then define objects resulting from a clicked x and y position on a spatial plot.

            Next, define a function to plot the spectra based on these two sets of x and y coordinates on the map. This will allow us to return spectra from a position we clicked on the image, and spectra where the mouse is currently hovering, allowing comparison of pixel reflectance values.

            In [24]:
            # Stream of X and Y positional data
            posxy = hv.streams.PointerXY(source=map, x=-61.833, y=-39.710) 
            clickxy = hv.streams.Tap(source=map, x=-61.833, y=-39.710) 
            
            # Function to build a new spectral plot based on mouse hover positional information retrieved from the RGB image using our full reflectance dataset
            def point_spectra(x,y):
                return ds_geo.sel(longitude=x,latitude=y,method='nearest').hvplot.line(y='reflectance',x='wavelengths',
                                                                                       color='#1b9e77', frame_width=400)
            # Function to build spectral plot of clicked location to show on hover stream plot
            def click_spectra(x,y):
                return ds_geo.sel(longitude=x,latitude=y,method='nearest').hvplot.line(y='reflectance',x='wavelengths',
                                                                                       color='#d95f02', frame_width=400)
            # Define the Dynamic Maps
            point_dmap = hv.DynamicMap(point_spectra, streams=[posxy])
            click_dmap = hv.DynamicMap(click_spectra, streams=[clickxy])
            
            # Plot the Map and Dynamic Map side by side
            (map + click_dmap*point_dmap)
            
            Out[24]:

            Contact Info:¶

            Email: LPDAAC@usgs.gov
            Voice: +1-866-573-3222
            Organization: Land Processes Distributed Active Archive Center (LP DAAC)¹
            Website: https://lpdaac.usgs.gov/
            Date last modified: 01-19-2022

            ¹Work performed under USGS contract G15PD00467 for NASA contract NNG14HH33I.